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Non uniform structures of the nucleon matter at subnuclear densities are numerically studied by 
means of the density functional theory with relativistic mean-fields coupled with the electric field. 

A particular role of the charge screening effects is demonstrated. 


I. INTRODUCTION 

There emerged many studies of the mixed phases at 
various first order phase transitions such as hadron-quark 
deconfinement transit icmjflS iiil , kaon condensa¬ 
tion Ej lajiiamE ll Jll4lTl5j[ T color superconductiv¬ 
ity liallalia. superfluidity in atomic traps [ill , nuclear 
pastaliflliJlH 0 0 111 EE El 0 111 etc. 

At very low densities, nuclei in matter are expected 
to form the Coulomb lattice embedded in the neutron- 
electron seas, that minimizes the Coulomb interaction 
energy. With an increase of the density, “nuclear pasta” 
structures emerge [20] : stable nuclear shape may change 
from droplet to rod, slab, tube, and to bubble. Pas¬ 
tas are eventually dissolved into uniform matter at a 
certain nucleon density below the saturation density, 
po — 0.16 fm -3 . Existence of pasta phases instead of 
the separated crystalline lattice of nuclei and the liq¬ 
uid npe phase would modify some important processes 
by changing the hydrodynamic properties and the neu¬ 
trino opacity in the supernova matter and in the matter 
of newly born neutron stars [^. Also the pasta phases 
may influence neutron star quakes and pulsar glitches via 
the change of mechanical properties of the crust matter 

El. 

A number of authors have investigated the low-density 
nuclear matter using various models [2(1 EH EE EE EE 
El EE El El EE Roughly speaking, the favorable nu¬ 
clear shape is determined by a balance between the sur¬ 
face and Coulomb energies. In most of the previous stud¬ 
ies the rearrangement effect of the density profile of the 
charged particles due to the Coulomb interaction is dis¬ 
carded. In Ref. the electron screening effect has been 
studied and it has been found that this effect is of minor 
importance. However, the rearrangement of the proton 
profiles as the consequence of the Coulomb repulsion was 
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not shown up in their model. 


A naive application of Gibbs conditions to separate 
bulk phases at the first order phase transitions, when 
one ignores the surface and Coulomb interaction, demon¬ 
strates a broad region of the structured mixed phase, cf. 
US However the charge screening effect (caused by 
the non-uniform charged particle distributions) should 
be very important when the typical structure size is of 
the order of the minimal Debye screening length in the 
problem. It may largely affect the stability condition of 
the geometrical structures in the mixed phases. We have 
been recently exploring the effect of the charge screening 
in the context of the various structured mixed phases 
S S S 0 El EE In fact, we have examined the 
mixed phase at the quark-hadron transition, kaon con¬ 
densation and of nuclear pasta, and found that in cases 
of the quark-hadron transition and kaon condensation 
the mixed phase might be largely limited by the charge 
screening and surface effects. 


Our purpose here is, following our preliminary study 
Ell, to investigate the nuclear pasta structures by 
means of a relativistic mean-field (RMF) model, which on 
the one hand does not need an introduction of the surface 
tension and on the other hand includes the Coulomb in¬ 
teraction in a proper way. We figure out how the charge 
screening effects modify the results obtained disregard¬ 
ing these effects. In Sec. EH we formulate the model 
and describe our numerical procedure. In Sec. Em we 
demonstrate the efficiency of the model in the descrip¬ 
tion of properties of finite nuclei. Then in Sec. II VI we de¬ 
scribe non-uniform pasta structures first at fixed proton 
to baryon number ratio that may have an application to 
the supernova matter and to the matter of a newly born 
hot protoneutron star. Then we investigate nuclear pasta 
at the beta-equilibrium, as they occur in cold neutron 
stars. In Sec. El we elucidate the effects of the surface 
and the charge screening. Then in Sec. EH we arrive at 
the conclusions. 
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II. DENSITY FUNCTIONAL THEORY WITH 
RELATIVISTIC MEAN FIELD 

A. Thermodynamic potential and equations of 
motion 


Following the idea of the density functional theory 
(DFT) with the RMF model §h|, we can formulate equa¬ 
tions of motion to study non-uniform nuclear matter nu¬ 
merically. The RMF model with fields of mesons and 
baryons introduced in a Lorentz-invariant way is rather 
simple for numerical calculations, but realistic enough to 
reproduce the bulk properties of nuclear matter. In our 
framework, the Coulomb interaction is properly included 
in equations of motion for nucleons, electrons, and me¬ 
son mean-fields, and we solve the Poisson equation for 
the Coulomb potential Vcoui self-consistently with them. 
Thus the baryon and electron density profiles, as well 
as the meson mean-fields, are determined in a way fully 
consistent with the Coulomb potential. 

Note that our framework can be easily extended to 
other situations; for example, if we take into account kaon 
either pion condensations, which are likely realized in 
a high-density region, we should only add the relevant 
meson field terms. In Ref. 0 we have included the 
kaon degree of freedom in such a treatment to discuss 
kaon condensation in high density regime. 

To begin with, we present the thermodynamic poten¬ 
tial for the system of neutrons, protons and electrons 
with chemical potentials, /i a (a = n,p,e), respectively; 

12 = 12^ T T H e , (1) 

where 

^ Jd 3 r 

a=p,n J 

with the local Fermi momenta, kp ta (r) (a = n,p), for 
nucleons, 


f>krp ( 
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for the scalar (a) and vector mean-fields (u>o,Ro) and 
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for electrons and the Coulomb potential, Vcoui(r), where 
i/ p (r) = p p + Vcoui(r) - g W Arwo(r) - g pN Ro(r), v n (r) = 
Pn - guNV o(r) + g pN R 0 (r), m* N ( r) = m N - g aN <j(r), 
and the nonlinear potential for the scalar field, U(cr) = 
\bv%N{gaNv) 3 + \c(g a N<r) A . Temperature is assumed to 
be zero in the present study. 


Here we use the local-density approximation for nucle¬ 
ons and electrons. Strictly speaking, the introduction of 
the density variable is meaningful, if the typical length 
of the nucleon density variation inside the structure is 
larger than the inter-nucleon distance. We must also 
bear in mind that for small structure sizes, quantum ef¬ 
fects become prominent which we disregarded. For the 
sake of simplicity we also discard nucleon and electron 
density derivative terms Is3| . In the case when we sup¬ 
press derivative terms of nucleon densities they follow 
changes of the other fields that have derivative terms. In 
our case these fields are the meson mean-fields and the 
Coulomb field. Here we consider large-size pasta struc¬ 
tures and simply discard the density variation effect, as a 
first-step calculation, while it can be easily incorporated 
in the quasi-classical manner by the derivative expansion 
within the density functional theory [34|. We also could 
use the fact that the resulting Debye screening lengths of 
electrons and protons characterizing the Coulomb field 
profile are typically much larger than those for all meson 
mean-fields. Then we could reduce contribution of the 
latter to the surface tension term. If the nucleon (neu¬ 
tron and proton) length scales were shorter than those of 
changes of the meson mean-fields, one could simplify the 
problem by dropping them and introducing instead a sur¬ 
face tension term. This simplified treatment is discussed 
in detail elsewhere l35i| . In this paper we avoid this sim¬ 
plification and solve the coupled-channel problem for the 
meson mean-fields and the Coulomb field numerically. 
Parameters of the RMF model are set to reproduce satu¬ 
ration properties of nuclear matter: the minimum energy 
per baryon —16.3 MeV at p = po = 0.153 fm~ 3 , the in¬ 
compressibility K(po) = 240 MeV, the effective nucleon 
mass m* N {po) = 0.78mjv; mjv = 938 MeV, and the sym¬ 
metry energy coefficient a sym = 32.5 MeV. Coupling 
constants and meson masses used in our calculation are 
listed in Table I. 

From the variational principle = 0 (d>i = 

a,R 0 ,u 0 , Vcoui) and = 0 (a = n,p, e), we get the 

coupled equations of motion for the mean-fields and the 
Coulomb potential, 


V 2 cr(r) = m 2 a(r) + — 
da 

-gvNip^ir) + P ( p s) { r)), (5) 

V 2 w 0 (r) = m 2 w 0 (r) - g uN {p P {r) + p n {r)), (6) 
V 2 I?o( r ) = m 2 p R 0 (r) - g pN (p P (r) - p n (r)), (7) 

V 2 Ucoui(r) = 47re 2 p ch (r), (8) 


with the scalar densities pi\ r) (a = n,p), and the charge 
density, p c h(r) = p p ( r) + p e (r). Equations of motion for 
fermions yield the standard relations between the densi¬ 
ties and chemical potentials, 




Pb = ^k 2 n (r) + m* N {r) 2 
+ 5 ^Arw 0 (r) - g pN R 0 (r), 


( 9 ) 
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TABLE I: Parameter set used in RMF in our calculation. 


9&N 

g^N 

g P N 

b 

c 

m CT [MeV] 

rriuj [MeV] 

m p [MeV] 

6.3935 

8.7207 

4.2696 

0.008659 

0.002421 

400 

783 

769 



+S^ArtA)(r) +g P Ari?o(r) - Vboui(r), (10) 

Pe{ r) = ~{p e - Vcoui(r)) 3 /37T 2 , (11) 

where we have assumed that the system is in chemi¬ 
cal equilibrium among nucleons and electrons and intro¬ 
duced the baryon-number chemical potential p b and the 
electron-number chemical potential p, e . Note that first, 
the Poisson equation for the Coulomb held © is highly 
nonlinear in Vcoui( r )> since p c h(r) in r.h.s. includes it in 
a complicated way. The Coulomb potential always en¬ 
ters equations through the gauge invariant combinations 
Pe hcoui(r) and p p -t- Vcoui(r). 

B. Numerical procedure 

To solve the above coupled equations numerically, we 
use the Wigner-Seitz cell approximation: the whole space 
is divided into equivalent cells with a geometry. The 
geometrical shape of the cell changes: sphere in three- 
dimensional (3D) calculation, cylinder in 2D and slab in 
ID, respectively. Each cell is globally charge-neutral and 
all the physical quantities in a cell are smoothly con¬ 
nected to those of the next cell with zero gradients at the 
boundary. Every point inside the cell is represented by 
the grid points (-/V gr id ss 100 ) and the differential equa¬ 
tions for fields are solved by the relaxation method for 
a given baryon-number density under the constraints of 
the global charge neutrality. 

To illustrate how to numerically solve equations of 
motion for the mean-fields, let us consider, for simplic¬ 
ity, two fields fi(r), / 2 (r) and their coupled Poisson-like 
equations under 3D calculation, 

V 2 /i(r) = mi 2 /i(r) + WM/ 1 ,/ 2 ], 

V 2 / 2 (r) = m 2 2 / 2 (r) + W 2 [/i,/ 2 ], (12) 

where Wi (i = 1 , 2 ) are functions of the fields /1 and / 2 . 
Introducing a relaxation “time” t artificially, we solve the 
equation, 

= d - m 2 fi(r;t) - Wj[/i,/ 2 ]) . 

(13) 

If the coefficients Ci are appropriately chosen, the above 
fi(r;t ) will converge to be constant in time and we get 
the solution of Eq. ©. 

The profiles of the nucleon densities are solved with the 
help of the “local chemical potentials” /i a (r) (a = n,p), 


being different from the constant chemical potentials 
which we have initially introduced. Assuming /z a (r) be¬ 
ing an increasing function of the neutron or proton num¬ 
ber density p a (r) in Eqs. m and m, the relaxation 
equation for the neutron or proton density profile, 

dPa Q t ,tS> = f ) Paid f)V 2 /z a (r; t), (14) 

is solved to get rid of the spatial dependence of the lo¬ 
cal chemical potentials p, a [r\t). The coefficients c Q (r;f) 
(a = n,p) are not constant so as to conserve the to¬ 
tal proton and neutron numbers. When we impose the 
beta-equilibrium condition, proton and neutron densities 
are adjusted to achieve p n (r) = p p {r) + p e {f)- Finally 
we get the density profiles p n (r) and p p (r ) relating to the 
constant chemical potentials p n (i") = p n and p p (r) = p p . 
Although the basic idea is to attain the constant chemi¬ 
cal potentials, p a (r) — p a {a = n,p) at the convergence, 
there is an exception: when there are some regions where 
Pair) = 0 , the local chemical potentials p a {r) are larger 
than the constant value in the regions where p Q (r) 7 ^ 0 . 

The electron density profile p e (r) is calculated directly 
from Eq. CD- The value of p e is adjusted at any time 
step to maintain the global charge neutrality: we de¬ 
crease p e when the total charge in a cell is positive and 
increase when it is negative. 

All the above relaxation procedures are performed si¬ 
multaneously. 

III. BULK PROPERTIES OF FINITE NUCLEI 

Before applying our framework to the problem of the 
pasta phase in nucleon matter, we check how it works 
to describe finite nuclei. In this calculation, for sim¬ 
plicity, we assume the spherical shape of nuclei. The 
electron density is set to be zero. Therefore neither the 
global charge neutrality condition nor the local charge- 
neutrality condition is imposed. 

In Fig. ^ (left panel) we show the density profiles of 
some typical nuclei. One can see how well our frame¬ 
work may reproduce the density profiles. To get a still 
better fit, especially around the surface region, we might 
need to include the derivative terms of the nucleon den¬ 
sities, as we have already remarked. Fine structures seen 
in the empirical density profiles, which may come from 
the shell effects (see, e.g., a proton density dip at the 
center of a light 16 O nucleus), cannot be reproduced by 
the mean-field theory. The effect is seen of the rearrange¬ 
ment of the proton density distribution in heavy nuclei. 
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FIG. 1: Properties of finite nuclei obtained with the present 
RMF model. Left: the density profiles of typical nuclei. The 
proton number densities (solid curves) are compared with the 
experiment. Right: binding energy per nucleon and the pro¬ 
ton number ratio of finite nuclei. 


Protons repel each other, which enhances their contam¬ 
ination near the surface of heavy nuclei. This effect is 
analogous to the charge screening effect for the Coulomb 
potential in a sense that the proton distribution is now 
changed not on the scale of the nuclear radius, but on 
another length scale, that we will call the proton De¬ 
bye screening length, see Eq. ifTTill below. It gives rise 
to important consequences for the pasta structures since 
typically the proton Debye screening length is less than 
the droplet size. The optimal value of the proton ( Z ) to 
the total baryon (A) number ratio Y p = Z/A is obtained 
by imposing the beta equilibrium condition for a given 
baryon number. Figure ^ (right panel) shows the baryon 
number dependencies of the binding energy per baryon 
and the proton number ratio. We can see that the bulk 
properties of finite nuclei (density, binding energy, and 
proton to baryon number ratio) are satisfactorily repro¬ 
duced for our present purpose. 


Note that in our framework we must use a sigma mass 
to ct = 400 MeV [36|, a slightly smaller value than that 
one usually uses, to get an appropriate fit. If we used a 
popular value rri„ ~ 500 MeV, finite nuclei would be 
over-bound by about 3 MeV/A The actual value of 
the sigma mass (as well as the omega and rho masses) 
has little relevance for the case of infinite nucleon mat¬ 
ter, since it enters the thermodynamic potential only in 
the combination C a = g a N/in a . However meson masses 
are important characteristics of finite nuclei and of other 
non-uniform nucleon systems, like those in pasta. The 
effective meson mass characterizes the typical scale for 
the spatial change of the meson field and consequently it 
affects the value of the effective surface tension |35j . 


IV. NON-UNIFORM STRUCTURES IN 
NUCLEON MATTER 


A. Nucleon matter at fixed proton number ratios 


First, we are concentrated on the discussion of the be¬ 
havior of the nucleon matter at a fixed value of the pro¬ 
ton number ratio Y p . Particularly, we explore the proton 
number ratios Y p = 0.1, 0.3, and 0.5. The cases Y p = 0.3 

0.5 should be relevant for the supernova matter and for 
newly born neutron stars. Figure [2 shows some typical 
density profiles inside the Wigner-Seitz cells. The geo¬ 
metrical dimension of the cell is denoted as “3D” (three- 
dimensional sphere), etc. The horizontal axis in each 
panel denotes the radial distance from the center of the 
cell. The cell boundary is indicated by the hatch. From 
the top to the bottom the configuration corresponds to 
droplet (3D), rod (2D), slab (ID), tube (2D), and bub¬ 
ble (3D). The nuclear “pasta” structures are clearly man¬ 
ifested. For the lowest Y p case ( Y p = 0.1), the neutron 
density is finite at any point: the space is filled by dripped 
neutrons; the neutron-drip value of Y p is, e.g., around 
0.26 in our 3D calculation. For a higher Y p , the neutron 
density drops to zero outside the nucleus. The proton 
number density always drops to zero outside the nucleus. 
We can see that the charge screening effects are pro¬ 
nounced. Due to the spatial rearrangement of electrons 
the electron density profile becomes no more uniform. 
This non-uniformity of the electron distribution is more 
pronounced for a higher Y p and a higher density. Protons 
repel each other. Thereby the proton density profile sub¬ 
stantially deviates from the step-function. The proton 
number is enhanced near the surface of the nucleus. 

The equation of state (EOS) for the sequence of ge¬ 
ometric structures is shown in Fig. 0 (top panels) as a 
function of the averaged baryon-number density. Note 
that the energy E — also includes the kinetic energy 
of electrons, which makes the total pressure positive. The 
lowest-energy configurations are selected among various 
geometrical structures. The most favorable configuration 
changes from the droplet to rod, slab, tube, bubble, and 
to the uniform one (the dotted thin curve) with an in¬ 
crease of density. The appearance of non-uniform struc¬ 
tures in matter results in a softening of EOS: the energy 
per baryon gets lower up to about 15 MeV/A compared 
to the uniform matter. 

The middle panels in Fig. [Hare partial pressures with¬ 
out electron contribution versus averaged baryon number 
density. If electron partial pressure is included, the total 
pressure becomes positive at all densities. 

The bottom panels in Fig.[o]are the cell radii R ce \\ and 
nuclear radii Rd versus averaged baryon number density. 
The radius Rd is defined by way of a density fluctuation 
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r [fm] r [fm] r [fm] 


FIG. 2: Examples of the density profiles in the cell for symmetric nuclear matter with Y P = 0.5 (left panel) and for asymmetric 
matter with Y v = 0.3 (center panel) and 0.1 (right panel). 



Pb [fm 3 ] ps [fm 3 ] Pb [fm 3 ] 


FIG. 3: Binding energy per nucleon, baryon partial pressure, and the cell and nuclear sizes for symmetric nuclear matter with 
Y p = 0.5 (left panel), and for asymmetric matter with Y' p = 0.3 (center panel) and 0.1 (right panel). 
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(for droplet, rod, and slab) 



(15) 

(for tube and bubble) 


where the bracket “()” indicates the average along the 
radial (for 3D and 2D cases) or perpendicular (ID) direc¬ 
tion in the cell. Dashed curves show the Debye screening 
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lengths of electron and proton calculated as 


A^ = ( —47re 


i d p\ 


dp e 


\(p) __ 
a d ~ 


47re' 


av \ —1/2 
2 '"t'p 


dPp 

dp p 


(16) 

(17) 


where is the proton number density averaged inside 
the nucleus (the region with finite p p ) and pf* is the 
electron charge density averaged inside the cell. Actu¬ 
ally doing more carefully we should introduce four De¬ 
bye screening lengths A/'^ and A^,’ > ' 1 with a separate 
averaging for the interior and the exterior of the nuclei. 
However we observe that the proton number density is al¬ 
ways zero in the exterior region and A^’^ = oo thereby. 

For electrons A^’ < ' 1 and A^’ > ' 1 are in general different 
but both being large and of the same order of magni¬ 
tude in the pasta case under consideration. Therefore 
we actually do not need a more detailed analysis of these 
quantities. Note that these values are obviously gauge in¬ 
variant. Numerically, the cell radii R ce n for droplet, rod, 
and slab configurations at Y p = 0.5 and 0.3 were proven 
to be close to the electron screening length. For the tube, 
f?ceii is larger than A}y. For Y p = 0.1, in all cases i? ce ii is 

(e) 

substantially smaller than X D and the electron screening 
should be much weaker, thereby. In all cases, except for 
bubbles (at Y p = 0.5 and 0.3), the structure radii Rd are 
smaller than X D . This means that the Debye screening 
effect of electrons inside these structures should not be 
pronounced. For bubbles at Y p = 0.5 and 0.3, X^f is 
substantially smaller than the cell size and the electron 
screening should be significant, see Fig. [5] below. For 
Y p = 0.5, 0.3, 0.1 in all cases (with the only exception 
Y p = 0.1 for slabs), the value X^’^ is shorter than Rd- 
Hence the density rearrangement of protons is essential 
for the pasta structures, as it is indeed seen from the Fig. 

13 

Knowing the baryon number density and the nuclear 
radius from Fig. [3] one may estimate the atomic number 
of the nucleus. In the case of droplets and for Y p = 0.5 the 
atomic number of the droplet is ~ 25 in the low density 
limit and ~ 65 at the maximum density of the droplet 
phase — 0-025 fm~ 3 . 


B. Nucleon matter in beta equilibrium 

Next, we consider the neutron star matter at zero tem¬ 
perature, and explore the non-uniform structures for the 
nucleon matter in beta equilibrium. Figure 0] shows the 
density profiles for different baryon number densities. 
The droplet structure itself is quite similar to the case 
of the fixed proton number ratio Y p = 0.1 considered 
above. The apparently different feature in this case is 
that only the droplet configuration appears as a non- 
uniform structure. It should be noticed, however, that 


the presence or absence of the concrete pasta structure 
sensitively depends on the choice of the effective interac¬ 
tion. 

In Fig. Owe plot the energy per baryon (top), the cell 
and nuclear sizes (middle), and the proton number ra¬ 
tio (bottom). The effect of the non-uniform structure on 
EOS (the difference between the energy of uniform mat¬ 
ter and that of non-uniform one) is small. However, the 
proton number ratio is significantly affected by the pres¬ 
ence of the pasta at lower densities. In the zero-density 
limit, the proton number ratio should converge to that of 
the normal nuclei. The droplet radius and the cell radius 
in the middle panel of Fig.0are always smaller than the 
electron Debye screening length X D . Thereby the effect 
of the electron charge screening is small. Unlike the fixed 
Y p cases, the droplet radius is comparable to the proton 
Debye screening length, which means that the effect of 
the proton rearrangement is not pronounced in this case. 
In fact, there is no enhancement of the proton number 
density near the surface in Fig.^J in contrast to Fig. [3 


V. COMPARISON WITH OTHER 
CALCULATIONS 

In this section we compare our DFT calculation with 
others to explore the effects of the surface, the charge 
rearrangement, and the fully consistent treatment of the 
density distribution. 

First let us focus on a very simplified treatment that 
has been used in the literature. We consider the bulk 
calculation supplemented by a simplified treatment of 
the finite-size effect. For description of the latter we in¬ 
troduce the surface tension and a bare (non screened) 
Coulomb interaction. This calculation assumes a sharp 
boundary between dense and dilute phases, uniform 
baryon density distribution inside each phase, and uni¬ 
form electron density distribution all over the cell. To 
further specify this approximation we use the term “no 
Coulomb + sharp surface”. We totally discard the 
Coulomb potential in equations of motion and drop the 
Poisson equation (“no Coulomb”) and we reduce the 
mean fields to their constant bulk values in the inte¬ 
rior and the exterior of the structure (“+ sharp sur¬ 
face”). The Coulomb energy, being evaluated with the 
step-function-like density profiles, and the surface energy, 
being expressed via the surface tension parameter r sur f 
are added to the total bulk energy. 

The volume fraction of each phase is simply calculated 
without taking into account of the finite-size effect (“bulk 
calculation”). Details of the “no Coulomb + sharp sur¬ 
face” calculation are presented in Appendix DTI 

Figure IH1 shows the EOS obtained by the “no Coulomb 
+ sharp surface” calculation performed at different val¬ 
ues of the surface tension. In this case for Y p = 0.5, 
the dilute phase includes no baryon. The value of the 
surface tension parameter r sur f ~ 1.03 MeV/fm 2 fits the 
liquid-drop binding energies of finite nuclei. Note that 
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FIG. 4: Density profiles in the cell for nuclear matter in 
beta equilibrium with baryon-number densities, 0.01, 0.03, 
0.05 and 0.07 fm~ 3 from top to bottom. 



FIG. 5: Binding energy (top), cell and nuclear sizes (mid¬ 
dle), and proton number ratio (bottom) in the cell for nuclear 
matter in beta equilibrium. 



FIG. 6: Bulk calculations with the surface tension parameter T su rf = 1.03 and 0.5 MeV/fm 2 . 


the appearance or disappearance of the pasta structure 
essentially depends on the value of the surface tension. 
With a larger value of the surface tension, the density re¬ 
gion of the pasta structure reduces and even some of the 
structures, e.g. “tube” and “bubble”, disappear. With 
a smaller value of the surface tension the region of the 
pasta structure broadens and all kinds of pasta struc¬ 
tures appear if we take r sur f < 0.3 MeV/fm 2 . However, 
if we put surface tension zero, the mixed phase reaches 
up from zero to the saturation density po without any 
specific geometry. Therefore, from the given example we 


see that the surface tension plays a crucial role in the ap¬ 
pearance of pasta structures. Remember that, in the case 
under consideration, the pasta structures are realized by 
a balance of the surface tension and the bare Coulomb 
interaction, which reads E sur f = 2i?coui, where E sur f is 
the surface energy and i?coui is the bare Coulomb en¬ 
ergy. Therefore, the Coulomb interaction is important, 
as well. Please also note that the surface tension intro¬ 
duced here simulates effects of the spatial changes of the 
meson mean-fields. In our “full calculation” the latter 
effects are taken into account explicitly whereas purely 



















r [fm] r [fm] r [fm] 


FIG. 7: Comparison of the density profiles for different treatments of the Coulomb interaction. From the left: full calculation, 
without electron screening, and “no Coulomb” calculation. The proton number ratio is Ip=0.5 for all cases. 



FIG. 8: Comparison of the density profiles for different treatments of the Coulomb interaction. From the left: “full” calculation, 
“no electron screening”, and “no Coulomb” calculation. The proton number ratio is Y p = 0.5 for all cases. 


“bulk calculations” completely disregard these effects. 

Next we compare three kinds of calculations with dif¬ 
ferent treatments of the Coulomb interaction. One is the 
“full calculation” which we have done here. The second is 
the calculation that disregards electron screening (“no e- 
screening”): a constraint is used that the electron density 
should be uniform. In this calculation the Coulomb po¬ 
tential Vcoui in Eq. ifTTl) is replaced by a constant Vo = 0, 

p e = ~{p e - Vo) 3 /3tt 2 . (18) 

In the full calculation the value of Vo is arbitrary, and one 


can take Vo for the sake of convenience, e.g. as Vo = 0, 
either set it equals the averaged value of Vc ou i (r) over 
the cell: recall that Vcoui either p e alone do not have a 
physical meaning but only the combination p e — Vcoui is 
meaningful due to the gauge invariance, cf. 4, l|. How¬ 
ever in the case “no e-screening” the gauge invariance 
is violated as it follows from Eq. m, since we replace 
Vcoui to Vo = 0 in the equation for the electron chem¬ 
ical potential but remain Vcoui in the equation for the 
proton chemical potential and thus in the expression for 
the proton number density. We do this procedure just 
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FIG. 9: Comparison of the phase diagrams between different 
treatments of the Coulomb interaction. 


to demonstrate the efficiency of the proton rearrange¬ 
ment artificially suppressing that of the electron one. The 
third calculation “no Coulomb” is performed by totally 
discarding the Coulomb potential Vboui in equations of 
motion. Accordingly, the Poisson equation is discarded 
as well. After getting the density profiles, the Coulomb 
energy, being evaluated using charge densities thus deter¬ 
mined, is added to the total energy. This calculation is 
similar to the above discussed “no Coulomb + sharp sur¬ 
face” calculation. The difference is that the effect of the 
density variation near the structure surface is automati¬ 
cally incorporated in explicit “no Coulomb” calculation, 
while in the above “no Coulomb + sharp surface” cal¬ 
culation this effect is hidden in the value of the surface 
tension. 

In Fig.0compared are the density profiles for different 
treatments of the Coulomb interaction. The left panel is 
the same as that in Fig. 0 It demonstrates the “full” 
calculation. It seems that there is almost no difference 
between the nucleon density of the “full” calculation and 
that of “no e-screening” calculation (center). The case 
of “no Coulomb” calculation (right), contrarily, shows 
a significant difference especially in the proton number 
density. The reason is simple: the electron Debye screen¬ 
ing length is large, whereas the proton Debye screening 
length is rather short. Thus the proton screening effects 
are much more pronounced than the electron ones. 

The EOS as a whole (upper panels in Fig. [5)l shows 
almost no dependence on the treatments of the Coulomb 
interaction. This agrees with a general statement that 
the variational functional is always less sensitive to the 
choice of the trial functions than the quantities linearly 
depending on them. Nevertheless, sizes of the cell and 
the nucleus (lower panels in Fig. |SJ) especially for tube 
and bubbles are different. In the cases of the “full cal¬ 
culation” and “no e-screening”, the cell radii of “tube” 
and “bubble” structures and that of “slab” structure get 
larger with an increase of density, while they are mono- 
tonically decreasing in the case of “no Coulomb” calcula¬ 
tion. We see almost no differences between the “full” and 


“no e-screening” calculations that again demonstrates 
relevance of the proton screening and weakness of the 
electron screening effects. The only significant difference 
remains for bubbles as seen from Fig. [U] The other ef¬ 
fect illustrated by Fig. |9] is a difference in the density 
range for each pasta structure. The “full” treatment of 
the Coulomb interaction slightly increases the region of 
the nuclear pasta. For Y p = 0.1 the differences between 
“full” and “no e-screening” calculations are completely 
washed out. 


VI. SUMMARY AND CONCLUDING 
REMARKS 

We have discussed the low-density nuclear matter 
structures, “nuclear pasta”, and elucidated the charge 
screening effect. Using a self-consistent framework based 
on density functional theory and relativistic mean fields, 
we took into account the Coulomb interaction in a proper 
way and numerically solved the coupled equations of mo¬ 
tion to extract the density profiles of nucleons. 

First we have checked how realistic our framework is by 
calculating the bulk properties of finite nuclei, as well as 
the saturation properties of nuclear matter, and found 
that it can describe both features satisfactorily. One 
could still improve the consideration fitting other exper¬ 
imental data. For example, we could more carefully fit 
different terms in the Weiczecker equation like the sur¬ 
face energy and the shell terms. For that we might be 
need an improvement of our relativistic mean field model 
that does not include the gradient proton and neutron 
density terms. 

In isospin-asymmetric nuclear matter for fixed proton 
to baryon number ratios, we have observed the “nuclear 
pasta” structures with various geometries at sub-nuclear 
densities. These cases are relevant for the discussion of 
the supernova explosions and for the description of the 
newly born neutron stars. The appearance of the pasta 
structures significantly lowers the energy, i.e. softens the 
equation of state, while the energy differences between 
various geometrical structures are rather small. The spa¬ 
tial rearrangement of the proton and electron charge den¬ 
sities (screening) affect the geometrical structures. 

By comparing different treatments of the Coulomb in¬ 
teraction, we have seen that the self-consistent inclusion 
of the Coulomb interaction changes the phase diagram. 
In particular the region of the pasta structure is broader 
for “full calculation” compared to that with simplified 
treatments of the Coulomb interaction which have been 
used in the previous studies. The effect of the rearrange¬ 
ment of the proton distributions on the structures is much 
more pronounced compared to the effect of the electron 
charge screening. The influence of the charge screening 
on the equation of state, on the other hand, was found 
to be small. 

We have also studied the structure of the nucleon mat¬ 
ter in the beta equilibrium. We have found that only one 
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type of structures is realized: proton-enriched droplets 
embedded in the neutron sea. No other geometrical 
structures like rod, slab, etc. appeared. 

In application to the newly formed neutron stars like 
in supernova explosions, finite temperature and neutrino 
trapping effects become important, as well as the dy¬ 
namics of the first order phase transition with formation 
of the structures. It would be interesting to extend our 
framework to include these effects. 
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APPENDIX A: “BULK” AND “NO COULOMB + 
SHARP SURFACE” CALCULATION OF LOW 
DENSITY NUCLEON MATTER 


into account the above conditions m, we obtain a set 

of pn , Pn l \ Pp\ Pp l \ = -P (II) , /, and the bulk 
energy density ebuik in each phase for given p n and p p . At 
this point ebuik does not include the surface and Coulomb 
contributions. If one cannot find the solution with finite 
Pn^ and pp 11 ), the proton or neutron density of the dilute 
phase is set to be zero. In this case the corresponding 
chemical potential is larger in the phase (II) than in the 
phase (I), and the complete set of the Gibbs conditions 
is not fulfilled. We leave off here the discussion of the 
“bulk calculation”. 

Now let us specify the “no Coulomb + sharp surface” 
calculation. To consider the structure of the mixed phase, 
the balance between the Coulomb interaction and the 
surface one should be taken into account. Introducing 
an adjusting parameter of the surface tension r SU rf, we 
calculate the surface energy density for the given geo¬ 
metrical dimension D: 


^surf 


Tsurf f D 

Rd 


(A2) 


where Rd is the droplet radius. The Coulomb energy 
density can be calculated [2(J as 


The bulk calculation proceeds like this (§, [jj (7, 16s 2fil : 
first consider two semi-infinite matters, (I) dense and (II) 
dilute phases, with a sharp boundary. The Coulomb and 
surface interactions are discarded for a while. Conditions 
of thermal equilibrium at zero temperature are imposed 
for pressure P and chemical potential p a (a — n,p ) be¬ 
tween two phases: 


$ 


2 - Df 1 - 2 / r> 


D-2 


P?> 

) 2 3 4 Rd 2 f 'R 
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D + 2 ’ 


By minimization of e sur f +£Coui in Rd (the relation e sur f = 
2 ecoui), we get 


p(b _ p( n )_ 
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All) 
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= „(II) 
r'n rn • 
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In Sec. II V Al the case is considered when the beta equi¬ 
librium is not imposed but the proton to baryon number 
ratio Y p is fixed. 

The averaged densities are p v = fp p ^ + (1 — f)p P ^ and 

Pn — fPn) + (1 — f)pn V> ■ Here / is the volume fraction 
of the phase (I). The chemical potentials are calculated 
using the RMF model presented in this paper. Taking 


Rd 


ecoui T e S urf 
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(A 6 ) 


Comparing the energy density of the uniform matter e 
and those of mixed phases ebuik + Curf + ecoui with differ¬ 
ent geometrical dimension D, we can determine the most 
favorable configuration and its energy density. 
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